# $ find . | grep "R$" | xargs grep -n "FeatureScatter" --color=auto
# seurat-4.1.0/R/visualization.R:1909:FeatureScatter <- function(



#' Scatter plot of single cell data
#'
#' Creates a scatter plot of two features (typically feature expression), across a
#' set of single cells. Cells are colored by their identity class. Pearson
#' correlation between the two features is displayed above the plot.
#'
#' @param object Seurat object
#' @param feature1 First feature to plot. Typically feature expression but can also
#' be metrics, PC scores, etc. - anything that can be retreived with FetchData
#' @param feature2 Second feature to plot.
#' @param cells Cells to include on the scatter plot.
#' @param shuffle Whether to randomly shuffle the order of points. This can be
#' useful for crowded plots if points of interest are being buried. (default is FALSE)
#' @param seed Sets the seed if randomly shuffling the order of points.
#' @param group.by Name of one or more metadata columns to group (color) cells by
#' (for example, orig.ident); pass 'ident' to group by identity class
#' @param cols Colors to use for identity class plotting.
#' @param pt.size Size of the points on the plot
#' @param shape.by Ignored for now
#' @param span Spline span in loess function call, if \code{NULL}, no spline added
#' @param smooth Smooth the graph (similar to smoothScatter)
#' @param slot Slot to pull data from, should be one of 'counts', 'data', or 'scale.data'
#' @param combine Combine plots into a single \code{\link[patchwork]{patchwork}ed}
#' @param plot.cor Display correlation in plot title
#' @param raster Convert points to raster format, default is \code{NULL}
#' which will automatically use raster if the number of points plotted is greater than
#' 100,000
#' @param raster.dpi Pixel resolution for rasterized plots, passed to geom_scattermore().
#' Default is c(512, 512).
#' @param jitter Jitter for easier visualization of crowded points
#'
#' @return A ggplot object
#'
#' @importFrom ggplot2 geom_smooth aes_string
#' @importFrom patchwork wrap_plots
#'
#' @export
#' @concept visualization
#'
#' @aliases GenePlot
#'
#' @examples
#' data("pbmc_small")
#' FeatureScatter(object = pbmc_small, feature1 = 'CD9', feature2 = 'CD3E')
#'
FeatureScatter <- function(
  object,
  feature1,
  feature2,
  cells = NULL,
  shuffle = FALSE,
  seed = 1,
  group.by = NULL,
  cols = NULL,
  pt.size = 1,
  shape.by = NULL,
  span = NULL,
  smooth = FALSE,
  combine = TRUE,
  slot = 'data',
  plot.cor = TRUE,
  raster = NULL,
  raster.dpi = c(512, 512),
  jitter = TRUE
) {
  #(A1) 默认使用全部细胞
  cells <- cells %||% colnames(x = object)
  
  #(A2) 如果 shuffle, 则打乱顺序
  if (isTRUE(x = shuffle)) {
    set.seed(seed = seed)
    cells <- sample(x = cells)
  }

  #(A3) meta.data 新增 ident 列，内容就是 sce@active.ident
  object[['ident']] <- Idents(object = object)

  #(A4) group.by 默认为 ident
  group.by <- group.by %||% 'ident'

  #(A5) 获取数据，共3列：2个基因，1个分组（也可能是多个分组: c("xx1", "xx2")）
  data <-  FetchData(
    object = object,
    vars = c(feature1, feature2, group.by),
    cells = cells,
    slot = slot
  )

  #(A6) 如果基因1在返回的数据框 列名没找到，报错
  # grepl() 返回逻辑值。对用户的输入要做各种检查，不能轻易相信。
  if (!grepl(pattern = feature1, x = colnames(x = data)[1])) {
    stop("Feature 1 (", feature1, ") not found.", call. = FALSE)
  }
  # 基因2同样
  if (!grepl(pattern = feature2, x = colnames(x = data)[2])) {
    stop("Feature 2 (", feature2, ") not found.", call. = FALSE)
  }

  # (A7) 转为数据框
  data <- as.data.frame(x = data)

  # (A8) 重新命名这2个基因名字，why? //todo
  # 难道基因还会被添加前缀吗？ 要再看 FetchData()，确实在源码解析 9-2(1)-A5: 
  # 如果指定的变量不在指定的 Assay 中，则在其他 assay 中找到时要加其Key前缀
  feature1 <-  colnames(x = data)[1]
  feature2 <-  colnames(x = data)[2]

  # (A9) 对 group.by 进行循环
  for (group in group.by) {
  	# 如果该列不是因子，变为因子
    if (!is.factor(x = data[, group])) {
      data[, group] <- factor(x = data[, group])
    }
  }

  # (A10) 对 group.by 进行循环
  # 如果是多个元素，则绘图结果也是多个组成的list
  plots <- lapply(
    X = group.by,
    FUN = function(x) {
      # 这个函数，是画图的主力函数
      SingleCorPlot(
        data = data[,c(feature1, feature2)],
        col.by = data[, x],
        cols = cols,
        pt.size = pt.size,
        smooth = smooth,
        legend.title = 'Identity',
        span = span,
        plot.cor = plot.cor,
        raster = raster,
        raster.dpi = raster.dpi,
        jitter = jitter
      )
    }
  )

  #(A11) 如果list长度为1，直接返回其内容
  if (isTRUE(x = length(x = plots) == 1)) {
    return(plots[[1]])
  }

  #(A12) 如果合并，则返回合并后的拼图
  if (isTRUE(x = combine)) {
    plots <- wrap_plots(plots, ncol = length(x = group.by))
  }

  return(plots)
}






# $ find . | grep "R$" | xargs grep -n "SingleCorPlot" --color=auto
# /seurat-4.1.0/R/visualization.R:6827:SingleCorPlot <- function(

例子: SingleCorPlot(FetchData(pbmc, vars=c("CD8A", "CD3D")))



globalVariables(names = '..density..', package = 'Seurat') # 这一行看不懂干啥的 //todo
#' A single correlation plot
#'
#' @param data A data frame with two columns to be plotted
#' @param col.by A vector or factor of values to color the plot by
#' @param cols An optional vector of colors to use
#' @param pt.size Point size for the plot
#' @param smooth Make a smoothed scatter plot
#' @param rows.highight A vector of rows to highlight (like cells.highlight in
#' \code{\link{SingleDimPlot}})
#' @param legend.title Optional legend title
#' @param raster Convert points to raster format, default is \code{NULL}
#' which will automatically use raster if the number of points plotted is
#' greater than 100,000
#' @param raster.dpi the pixel resolution for rastered plots, passed to geom_scattermore().
#' Default is c(512, 512)
#' @param plot.cor ...
#' @param jitter Jitter for easier visualization of crowded points
#'
#' @return A ggplot2 object
#'
#' @importFrom stats cor
#' @importFrom cowplot theme_cowplot
#' @importFrom RColorBrewer brewer.pal.info
#' @importFrom ggplot2 ggplot aes_string geom_point labs scale_color_brewer
#' scale_color_manual guides stat_density2d aes scale_fill_continuous
#' @importFrom scattermore geom_scattermore
#'
#' @keywords internal
#'
#' @export
#' 
SingleCorPlot <- function( #FeatureScatter() 的主力函数
  data,
  col.by = NULL,
  cols = NULL,
  pt.size = NULL,
  smooth = FALSE,
  rows.highlight = NULL,
  legend.title = NULL,
  na.value = 'grey50',
  span = NULL,
  raster = NULL,
  raster.dpi = NULL,
  plot.cor = TRUE,
  jitter = TRUE
) {
  #(A1)点的大小
  pt.size <- pt.size %||% AutoPointSize(data = data, raster = raster)

  #(A2) 是否栅格化
  # 如果超过10万个点，且 raster 不是F，则设置为允许T
  if ((nrow(x = data) > 1e5) & !isFALSE(raster)){
    message("Rasterizing points since number of points exceeds 100,000.",
            "\nTo disable this behavior set `raster=FALSE`")
  }
  raster <- raster %||% (nrow(x = data) > 1e5)

  #(A3) 如果设置了分辨率 
  if (!is.null(x = raster.dpi)) {
  	# 如果不是数字，或 长度不是2，则报错
    if (!is.numeric(x = raster.dpi) || length(x = raster.dpi) != 2)
      stop("'raster.dpi' must be a two-length numeric vector")
  }

  # (A4) 记录df的列名
  orig.names <- colnames(x = data)
  
  # (A5) 替换列名，记录在 新的列名中、新变量中
  # 替换列名中的-为.
  names.plot <- colnames(x = data) <- gsub(
    pattern = '-',
    replacement = '.',
    x = colnames(x = data),
    fixed = TRUE
  )
  # 替换':' 为 '.'，
  names.plot <- colnames(x = data) <- gsub(
    pattern = ':',
    replacement = '.',
    x = colnames(x = data),
    fixed = TRUE
  )
  
  # (A6) 如果少于2列，报错
  if (ncol(x = data) < 2) {
    msg <- "Too few variables passed"
    if (ncol(x = data) == 1) {
      msg <- paste0(msg, ', only have ', colnames(x = data)[1])
    }
    stop(msg, call. = FALSE)
  }

  # (A7) 如果需要相关系数，就计算；否则空字符串
  plot.cor <- if (isTRUE(x = plot.cor)) {
    round(x = cor(x = data[, 1], y = data[, 2]), digits = 2)
  }
  else(
    ""
  )

  # (A8) 如果设置了 rows.highlight
  # 粗看一遍，看完画图回头再看二遍：看每个参数的影响范围
  # 这次没有走这个if，下次遇到了再说第二遍
  if (!is.null(x = rows.highlight)) {
  	# 设置高亮，返回值是一个list。见下文。
    highlight.info <- SetHighlight(
      cells.highlight = rows.highlight, # 要高亮显示的行 cell id
      cells.all = rownames(x = data), #全部cid

      sizes.highlight = pt.size, #要高亮的尺寸
      cols.highlight = 'red', #要高亮的颜色

      col.base = 'black', #背景颜色
      pt.size = pt.size #背景点大小
    ) #"plot.order" "highlight"  "size"       "color" 
    # 背景色
    cols <- highlight.info$color
    # 染色顺序
    col.by <- factor(
      x = highlight.info$highlight,
      levels = rev(x = highlight.info$plot.order)
    )
    # 排序
    plot.order <- order(col.by)
    data <- data[plot.order, ]
    col.by <- col.by[plot.order]
  }

  # 如果 col.by 非空，
  if (!is.null(x = col.by)) {
    data$colors <- col.by
  }

  ######################
  # (B1)开始画图
  plot <- ggplot(
    data = data,
    mapping = aes_string(x = names.plot[1], y = names.plot[2])
  ) +
    labs(
      x = orig.names[1],
      y = orig.names[2],
      title = plot.cor,
      color = legend.title #这个值有作用吗？我测试没作用
    )

  # (B2)如果光滑化，使用 stat_density2d() 画图
  if (smooth) {
    # density <- kde2d(x = data[, names.plot[1]], y = data[, names.plot[2]], h = Bandwidth(data = data[, names.plot]), n = 200)
    # density <- data.frame(
    #   expand.grid(
    #     x = density$x,
    #     y = density$y
    #   ),
    #   density = as.vector(x = density$z)
    # )
    plot <- plot + stat_density2d(
      mapping = aes(fill = ..density.. ^ 0.25),
      geom = 'tile',
      contour = FALSE,
      n = 200,
      h = Bandwidth(data = data[, names.plot])
    ) +
      # geom_tile(
      #   mapping = aes_string(
      #     x = 'x',
      #     y = 'y',
      #     fill = 'density'
      #   ),
      #   data = density
      # ) +
      scale_fill_continuous(low = 'white', high = 'dodgerblue4') +
      guides(fill = FALSE)
  }

  # (B3)位置
  position <- NULL
  # 判断是否抖动
  if (jitter) {
    position <- 'jitter'
  } else {
    position <- 'identity'
  }
  
  # (B4)是否指定颜色
  # (C1)如果指定颜色了，这个颜色用到哪里了？就是 mapping 中的 data$colors列，前面添加的。
  if (!is.null(x = col.by)) {
  	# 栅格化
    if (raster) {
      plot <- plot + geom_scattermore(
        mapping = aes_string(color = 'colors'),
        position = position,
        pointsize = pt.size, #这一行参数名不同
        pixels = raster.dpi #多了这一行
      )
    # 如果不栅格化
    } else {
      plot <- plot + geom_point(
        mapping = aes_string(color = 'colors'),
        position = position,
        size = pt.size
      )
    }
  # (C2)如果没有指定颜色
  } else {
    if (raster) {
      plot <- plot + geom_scattermore(position = position, pointsize = pt.size, pixels = raster.dpi)
    } else {
      plot <- plot + geom_point(position = position, size = pt.size)
    }
  }

  # (B5) 点的颜色，数量要和 levels 一样多
  # 如果有 cols 参数
  if (!is.null(x = cols)) {
  	# 如果指定的cols参数是1个，且在 rownames(RColorBrewer::brewer.pal.info) 中
    cols.scale <- if (length(x = cols) == 1 && cols %in% rownames(x = brewer.pal.info)) {
      scale_color_brewer(palette = cols) #则定义为调色板
    } else {
      # 否则，直接使用这些颜色值，指定 na.value 的颜色
      scale_color_manual(values = cols, na.value = na.value)
    }

    # 增加颜色层
    plot <- plot + cols.scale

    # 如果高亮非空，则不显示颜色图例
    if (!is.null(x = rows.highlight)) {
      plot <- plot + guides(color = FALSE)
    }
  }

  # (B6) 添加主题，标题居中
  plot <- plot + theme_cowplot() + theme(plot.title = element_text(hjust = 0.5))

  # (B7) 如果有 span 参数
  if (!is.null(x = span)) {
  	# 则增加 LOESS 曲线
    plot <- plot + geom_smooth(
      mapping = aes_string(x = names.plot[1], y = names.plot[2]),
      method = 'loess',
      span = span
    )
  }


  return(plot)
}






(2) AutoPointSize()
# seurat-4.1.0/R/visualization.R:4038:AutoPointSize <- function(data, raster = NULL) {

# 自动获取点的大小。
# 如果栅格化，取1；否则，取 1583/nrow,1 中较小的。还是不懂为什么这样设置。
AutoPointSize <- function(data, raster = NULL) {
  return(ifelse(
    test = isTRUE(x = raster),
    yes = 1,
    no = min(1583 / nrow(x = data), 1)
  ))
}








(3) SetHighlight()
# /seurat-4.1.0/R/visualization.R:6716:SetHighlight <- function(

# 测试：
> table(pbmc$seurat_clusters)
  0   1   2   3   4   5   6   7   8 
711 480 472 344 279 162 144  32  14

rs=Seurat:::SetHighlight(
  cells.highlight = WhichCells(pbmc, expression= seurat_clusters == 3),
  cells.all = colnames(pbmc),
  
  sizes.highlight = 2,
  cols.highlight = "red",
  
  col.base="black",
  pt.size=1
)

str(rs) #返回值是一个list
#List of 4
# $ plot.order: chr [1:2] "Group_1" "Unselected"
# $ highlight : Factor w/ 2 levels "Group_1","Unselected": 2 1 2 2 2 2 2 2 2 2 ...
# $ size      : num [1:2638] 1 2 1 1 1 1 1 1 1 1 ...
# $ color     : chr [1:2] "black" "red"

names(rs) #[1] "plot.order" "highlight"  "size"       "color" 
table(rs$highlight)
# Group_1 Unselected 
#     344       2294






源码：
# Set highlight information
#
# @param cells.highlight Cells to highlight
# @param cells.all A character vector of all cell names
# @param sizes.highlight Sizes of cells to highlight
# @param cols.highlight Colors to highlight cells as
# @param col.base Base color to use for unselected cells
# @param pt.size Size of unselected cells
#
# @return A list will cell highlight information
# \describe{
#   \item{plot.order}{An order to plot cells in}
#   \item{highlight}{A vector giving group information for each cell}
#   \item{size}{A vector giving size information for each cell}
#   \item{color}{Colors for highlighting in the order of plot.order}
# }
#
SetHighlight <- function(
  cells.highlight,
  cells.all,

  sizes.highlight,
  cols.highlight,

  col.base = 'black',
  pt.size = 1
) {
  #(A1) 如果要高亮的细胞是 字符串
  if (is.character(x = cells.highlight)) {
  	# 封到 list 中
    cells.highlight <- list(cells.highlight)
  # 如果是 df，或者不是list
  } else if (is.data.frame(x = cells.highlight) || !is.list(x = cells.highlight)) {
  	# 转为 list
    cells.highlight <- as.list(x = cells.highlight)
  }

  #(A2) 对每个要高亮的细胞，得到的是一个list
  cells.highlight <- lapply(
    X = cells.highlight,
    FUN = function(cells) {
      # 如果是字符串
      cells.return <- if (is.character(x = cells)) {
      	# 就取在所有细胞内的部分
        cells[cells %in% cells.all]
      # 如果不是字符串
      } else {
      	# 转成数字，当下标使用
        cells <- as.numeric(x = cells)
        # 下标超过总细胞数的不要
        cells <- cells[cells <= length(x = cells.all)]
        # 取总细胞的这些下标
        cells.all[cells]
      }
      # 返回要高亮的cell id
      return(cells.return)
    }
  )

  #(A3) 过滤，去掉长度为0的list键值对
  cells.highlight <- Filter(f = length, x = cells.highlight)


  #(A4) 高亮的细胞的name
  # 如果list没有名字，则名字为 Group_数字下标，长度和高亮细胞等长
  names.highlight <- if (is.null(x = names(x = cells.highlight))) {
    paste0('Group_', 1L:length(x = cells.highlight))
  # 如果有名字，则使用其名字
  } else {
    names(x = cells.highlight)
  }

  #(A5) 高亮 尺寸
  # 生成等长度的 尺寸 向量
  sizes.highlight <- rep_len(
    x = sizes.highlight,
    length.out = length(x = cells.highlight)
  )

  #(A6) 高亮 颜色
  # 第一个元素是 col.base，其余是 等长度的 高亮颜色 向量
  cols.highlight <- c(
    col.base,
    rep_len(x = cols.highlight, length.out = length(x = cells.highlight))
  )

  #(A7) 点的大小，和所有细胞等长
  size <- rep_len(x = pt.size, length.out = length(x = cells.all))

  #(A8) 高亮， 是一个NA字符串，和所有细胞等长
  highlight <- rep_len(x = NA_character_, length.out = length(x = cells.all))

  #(A9) 如果 高亮细胞list 大于0
  if (length(x = cells.highlight) > 0) {
  	# 对每个高亮细胞list内的元素进行循环
    for (i in 1:length(x = cells.highlight)) {
      # 第一个 list 中的 cell id
      cells.check <- cells.highlight[[i]]
      # 查他们在全部cell中的下标
      index.check <- match(x = cells.check, cells.all)

      # 该下标位置的 高亮 值，等于 高亮细胞的名字
      highlight[index.check] <- names.highlight[i]
      # 该下标位置的 尺寸 值，等于 高亮细胞的尺寸
      size[index.check] <- sizes.highlight[i]
    }
  }

  #(A10) 排序，默认是升序；设置 NA 放到最后
  plot.order <- sort(x = unique(x = highlight), na.last = TRUE)
  # 把NA值替换为 'Unselected'
  plot.order[is.na(x = plot.order)] <- 'Unselected'

  # (A11) 替换NA为 'Unselected'
  highlight[is.na(x = highlight)] <- 'Unselected'

  # 变为因子
  highlight <- factor(x = highlight, levels = plot.order)

  return(list(
    plot.order = plot.order,
    highlight = highlight,
    size = size,
    color = cols.highlight
  ))
}












# Bandwidth
# seurat-4.1.0/R/visualization.R:5320:Bandwidth <- function(data) {

测试：
> Seurat:::Bandwidth(FetchData(pbmc, vars=c("CD3D", "CD8A")))
[1] 0.5838914


# Calculate bandwidth for use in ggplot2-based smooth scatter plots
#
# Inspired by MASS::bandwidth.nrd and graphics:::.smoothScatterCalcDensity
#
# @param data A two-column data frame with X and Y coordinates for a plot
#
# @return The calculated bandwidth
#
#' @importFrom stats quantile var
#
Bandwidth <- function(data) {
  # diff返回下一行减去上一行的差。
  r <- diff(x = apply(
    X = data,
    MARGIN = 2,
    FUN = quantile, #求每一列的分位数，去掉na，不要名字
    probs = c(0.05, 0.95),
    na.rm = TRUE,
    names = FALSE
  ))
  # 两个的差的绝对值，除以 0.34
  h <- abs(x = r[2L] - r[1L]) / 1.34
  # 这个差如果等于0，则赋值为1
  h <- ifelse(test = h == 0, yes = 1, no = h)
  
  # 计算返回值，这几个数怎么定的呢? 黑魔法！
  bandwidth <- 4 * 1.06 *
    min(sqrt(x = apply(X = data, MARGIN = 2, FUN = var)), h) * # 最小值: 每列方差，前面的h
    nrow(x = data) ^ (-0.2)

  return(bandwidth)
}


